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We calculate the next-to-next-to-leading-order 0(ai) one-loop squared corrections to the produc- 
tion of heavy-quark pairs in the gluon-gluon fusion process. Together with the previously derived 
results on the qq production channel, the results of this paper complete the calculation of the one- 
loop squared contributions of the next-to-next-to-leading-order O(a^) radiative QCD corrections 
to the hadroproduction of heavy flavors. Our results, with the full mass dependence retained, are 
presented in a closed and very compact form, in dimensional regularization. 
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I. INTRODUCTION 

It has been already 20 years since the next-to-leading- 
order (NLO) corrections to the hadroproduction of heavy 
flavors were first presented in the seminal work [l[ . These 
results were confirmed yet in another seminal work 

In the past few years there was much progress in de- 
scribing the experimental results on heavy-flavor produc- 
tion. For instance, in a recent work [3| it was shown 
that a NLO analysis of the transverse-momentum distri- 
butions does in fact properly describe the latest bottom 
quark production data [J] in a surprisingly large kinemat- 
ical range. The improvement in the theoretical prediction 
is mainly due to advances in the analysis of parton dis- 
tribution functions and the QCD coupling constant. We 
also point out the progress in dealing with numerically 
large mass logarithms that spoil the convergence of the 
perturbative expansion in the high energy (or small mass) 
asymptotic domain. In this respect we mention the work 
Q where also charm pair production is reconciled with 
experimental data. Data on top-quark pair production 
also agrees with the NLO prediction within theoretical 
and experimental errors (see e.g. Ref. (6(). However, 
in all of these NLO calculations there remains, among 
others, the problem that the renormalization and factor- 
ization scale dependences render the theoretical predic- 
tions to have much larger uncertanties than today's stan- 
dards require. This calls for a next-to-next-to-leading- 
order (NNLO) calculation of heavy-quark production in 
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hadronic collisions. In fact, the scale dependence of the 
theoretical prediction is expected to be considerably re- 
duced when NNLO partonic amplitudes are folded with 
the available NNLO parton distributions. For example, 
by approximating the NNLO corrections with the fixed- 
order expansion of the next-to-leading-log prediction, one 
finds a projected NNLO scale uncertainty of about 3% 0] , 
which is below the parton distribution uncertainty, and 
in line with the anticipated experimental error. 

Recently there was much activity in the phenomenol- 
ogy of hadronic heavy-quark pair production in connec- 
tion with the Tevatron and the CERN Large Hadron Col- 
lider (LHC), which had its start-up this year. There will 
be much experimental effort dedicated to the discovery 
of the Higgs boson. There will also be studies of the 
copious production of top quarks and other heavy parti- 
cles, which serve as a background to Higgs boson searches 
as well as to possible new physics beyond the standard 
model. Therefore, it is mandatory to reduce the the- 
oretical uncertainty in phcnomcnological calculations of 
heavy-quark production processes as much as possible. 

Several years ago the NNLO contributions to hadron 
production were calculated by several groups in massless 
QCD (see e.g. Ref. Q and references therein). The com- 
pletion of a similar program for processes that involve 
massive quarks requires much more dedication, since the 
inclusion of an additional mass scale dramatically com- 
plicates the whole calculation. 

At the lower energies of Tevatron II, top-quark pair 
production is dominated by qq annihilation (85%). The 
remaining 15% comes from gluon fusion. At the higher 
energies of the LHC, gluon fusion dominates the produc- 
tion process (90%) leaving 10% for qq annihilation (per- 
centage figures from Ref. Q). This shows that both qq 
annihilation and gluon fusion have to be accounted for in 
the calculation of top-quark pair production. Since gluon 
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FIG. 1: Exemplary gluon fusion diagrams for the NNLO calculation of heavy-hadron production. 



fusion makes up the largest part of the heavy-quark pair 
production cross section at the LHC it is important to 
reduce renormalization and factorization scale uncertain- 
ties in the gluon fusion process as much as possible in 
view of the fact that the large uncertainties in the glu- 
onic parton distribution functions translate to large cross 
section uncertainties at the LHC. 

There are four classes of contributions that need to 
be calculated for the NNLO corrections to the hadronic 
production of heavy-quark pairs. In Fig. [1] we show one 
generic diagram each for the four classes of contribu- 
tions that need to be calculated for the NNLO correc- 
tions to the gluon-initiated hadroproduction of heavy fla- 
vors. The first class involves the pure two-loop contribu- 
tion QIa)], which has to be folded with the leading-order 
(LO) Born term. The second class of diagrams [TJb)] 
consists of the so-called one-loop squared contributions 
(also called loop-by-loop contributions) arising from the 
product of one-loop virtual matrix elements. This is the 
topic of the present paper. Further, there are the one- 
loop gluon emission contributions QIc)] that are folded 
with the one-gluon emission graphs. Finally, there are 
the squared two-gluon emission contributions [2(d)] that 
are purely of tree type. The corresponding graphs for the 
quark-initiated processes are not displayed. 

Bits and pieces of the NNLO calculation for hadropro- 
duction of heavy flavors are now being assembled. In 
this context we would like to mention the recent two- 
loop calculation of the heavy-quark vertex form factor @ 
that can be used as one of the many building blocks in 
the first class of processes. There is also a very promising 
numerical approach applied to the calculation of the pure 
two- loop diagrams [l0( • Recently, an analytic calculation 
of a subclass of the two-loop contributions to qq — > QQ 
was published [HI- The authors of Ref. [3 have cal- 
culated the NLO corrections to it+jet production with 
contributions from the third class of diagrams. However, 
this result needs further subtraction terms in order to 
allow for an integration over the full phase space. We 
would also like to mention the recent work on the two- 
loop virtual amplitudes that are valid in the domain of 
high energy asymptotics, where the heavy-quark mass is 
small compared to the other large scales. In this calcu- 
lation [13j ]. mass power corrections are left out, and only 



large mass logarithms and finite terms associated with 
them are retained. Much work was also done in relation 
to the resummation of soft contributions. In this respect 
we refer the reader to recent publications where some 
different approaches to the resummation are advocated 

The authors of the present paper have been involved in 
a systematic effort to calculate all the contributions from 
the second class of processes, i.e. the one-loop squared 
contributions. The NNLO one-loop squared amplitudes 
for the quark-initiated process were recently presented in 
Ref. In this paper, we report on a calculation of the 
NNLO one-loop squared matrix elements for the process 
99 ~* QQ- The calculation is carried out in dimensional 
regularization [l6| with space-time dimension n — 4 — 2e. 
We mention that we have presented closed-form, one- 
loop squared results for heavy-quark production in the 
fusion of real photons in Ref. [13]. With the present 
paper the program of calculating the one-loop squared 
contributions to heavy-quark pair hadroproduction has 
now been completed. 

Let us briefly describe some of the main features of the 
calculation of the one-loop squared contributions. The 
highest singularity in the one-loop amplitudes arises from 
infrared (IR) and mass singularities (M) and is thus, in 
general, proportional to (1/e 2 ). This in turn implies that 
the Laurent series expansion of the one-loop amplitudes 
has to be taken up to 0(e 2 ) when calculating the one- loop 
squared contributions. In fact, it is the 0(e 2 ) terms in the 
Laurent series expansion that really complicate things 
[HI], since the 0(s 2 ) contributions in the one-loop am- 
plitudes involve a multitude of multiple polylogarithms 
of maximal weight and depth 4 [ll|. All scalar master 
integrals needed in this calculation have been assembled 
in Refs. [H, Eil- Reference [H| gives the results in terms 
of so-called L functions, which can be written as one- 
dimensional integral representations involving products 
of log and dilog functions, while Ref. [n| gives the re- 
sults in terms of multiple polylogarithms. The divergent 
and finite terms of the one-loop amplitude for gg — ► QQ 
were given in Ref. 20]. The remaining O(s) and 0(e 2 ) 
amplitudes have been written down in Ref. (21| . We shall 
rewrite these matrix elements in a representation more 
suitable for the purposes of the present application. 
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FIG. 2: The i-, M-, and s-channel LO graphs contributing to the gluon (curly lines) fusion amplitude. The thick solid lines 
correspond to the heavy quarks. 



In our presentation, we shall make use of our nota- 
tion for the coefficient functions of the relevant scalar 
one-loop master integrals calculated up to 0{e 2 ) in 
Refs. [IH For the case of gluon-gluon and quark- 

antiquark collisions, one needs all the scalar integrals 
derived in Refs. [H, [l9j], e.g. the one scalar one-point 
function A, the five scalar two-point functions B\, B2, 
B3, -B4, and B5, the six scalar three-point functions 
C\, C2, C3, C4, C5, and 6*6, and three scalar four-point 
functions D\,D2, and D3. Taking the complex scalar 
four-point function D2 as an example, we define succes- 
sive coefficient functions D 2 ' for the Laurent series ex- 
pansion of Z?2- One has 

D 2 = iC£ (m 2 ){l^- 2 ) + Id™ +sDP 

+e 2 D { 2 2) +0(e 3 )}, (1.1) 

where C £ (m 2 ) is defined by 

We use this notation for both the real and imaginary 
parts of D2, i.e. for ReZ?2 and Im_D2- Similar expansions 
hold for the scalar one-point function A, the scalar two- 
point functions Bi, the scalar three-point functions Ci, 
and the remaining four-point functions Di. The coeffi- 
cient functions of the various Laurent series expansions 
were given in Ref. [18| in the form of so-called L func- 
tions, and in Ref. in terms of multiple poly logarithms 
of maximal weight and depth 4. It is then a matter of 
choice which of the two representations are used for the 
numerical evaluation. The numerical evaluation of the L 
functions in terms of their one-dimensional integral rep- 
resentations is quite straightforward using conventional 
integration routines, while there exists a very efficient al- 
gor ithm to numerically evaluate multiple polylogarithms 

Let us briefly summarize the main features of 
the scalar master integrals. The master inte- 
grals A, B\, B 3 , B 4} C2, C3, and D 3 are real, whereas 
B 2 , B 5 , Ci, C 4 , C 5 , C 6 , Di, and D 2 are complex. From 
the form (AB* + BA*) = 2(ReAReB + JmAImB) it is 
clear that the imaginary parts of the master integrals 



must be taken into account in the one-loop squared con- 
tribution. The master integrals B2, -B5, C\, C4, C5, and 
Cq are (t «-> u) symmetric, where the kinematic variables 
t and u are defined in Sec. HH 

This paper is organized as follows. Section. |TT] con- 
tains an outline of our general approach and discusses 
renormalization procedures. Section. IHII presents LO and 
NLO results for the gluon fusion subprocess. In Sec. IIVI 
one finds a discussion of the singularity structure of the 
NNLO squared matrix element for the gluon fusion sub- 
process. In Sec. |V]we discuss the structure of the finite 
part of our result. Our results are summarized in Sec. lVIl 
In the Appendices, we present expressions for various co- 
efficients that are used in Sec.|nT]to write down the NLO 
result. 

II. NOTATION AND RENORMALIZATION 

Heavy-flavor hadroproduction proceeds through two 
partonic subprocesses: gluon fusion and light-quark- 
antiquark annihilation. The first subprocess is the most 
challenging one in QCD from a technical point of view. 
It has three production topologies already at the Born 
level (see Fig. [2]). The second subprocess, where there is 
only one topology at the Born level, was considered in 
Ref. Irrespective of the partons involved, the gen- 

eral kinematics is, of course, the same in both processes. 
In particular, for gluon fusion, Fig. [2l we have 

g(pi)+g(p2)^Q(P3) + Q(P4), (2.1) 

The momentum flow directions correspond to the phys- 
ical configuration, e.g. pi and P2 are ingoing whereas p% 
and P4 are outgoing. With m being the heavy-quark 
mass, we define 

s=( Pl + P 2) 2 , t = T -m 2 = (pi-p 3 ) 2 -to 2 , 

u=U -m 2 = (p 2 -p 3 ) 2 -m 2 , (2.2) 

so that one has the energy-momentum conservation rela- 
tion s + t + u = 0. 

We also introduce the overall factor 

C = (giC E (m 2 )) 2 , (2.3) 

where g s is the renormalized strong-coupling constant 
and C £ (m 2 ) is defined in Eq. (|1.2|) . 
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As was shown e.g. in Refs. [2CJ, |21| the self-energy and 
vertex diagrams contain ultraviolet (UV), infrared and 
collinear (IR/M) poles after heavy-mass renormalization. 
The UV poles need to be regularized. 

Our renormalization procedure is carried out in a 
mixed renormalization scheme. When dealing with mass- 
less quarks, we work in the modified minimal-subtraction 
(MS) scheme, while heavy quarks are renormalized in 
the on-shell scheme defined by the following conditions 
for the renormalized external heavy-quark self-energy 
graphs: 



Sr(^)|jC=m = 0, 



d_ 



E r(i>)|j/= 



0. 



(2.4) 



In the on-shell scheme, the first condition in Eq. 
ensures that the heavy-quark mass is the pole mass. 

For completeness, we list the set of one-loop renormal- 
ization constants used in this paper. One has 



Z m = l-g 2 C F C £ (m 2 ) 
Z2 = Z m , 



3- 2e 
e(l-2e) ! 



(2.5) 



Zip = Z 2 - »±N c C £ (p 2 ), 

£ 

Z lf = l-^N C C £ (^ 2 ), 

Z 3 = l + & {(~N C - ^)C £ ( M 2 ) - 2 -C e {m 2 ) 
= 1 + - 2N c )C e {n 2 ) - lc £ (m 2 ) 

£ I 3 



Z g = l- 9 l\^C e {f)-\c e {m 2 ))> 

with Pq = (UNc — 2ni)/3 being the first coefficient of 
the QCD beta function, ni the number of light quarks, 
Cf = 4/3, and Nc — 3 the number of colors. The arbi- 
trary mass scale /i is the scale at which the renormaliza- 
tion is carried out. The above renormalization constants 
renormalize the following quantities: Z\ for the three- 
gluon vertex, Z m for the heavy-quark mass, Z 2 for the 
heavy-quark wave function, Zip for the {QQg) vertex, 
Zi f for the (qqg) vertex, Z3 for the gluon wave func- 
tion and Z g for the strong-coupling constant a s . For 
the massless quarks, there is no mass and wave function 
renormalization. 

Let us sketch the two alternative ways of getting the 
final one-loop-renormalized amplitude from the mass- 
renormalized amplitude: 

i) Take the given mass-renormalized matrix element or 
the square of that matrix element and multiply all the 
self-energy graphs by a factor 1/2. Then renormalize the 
coupling constant in the LO Born amplitude. 

ii) Take the given mass-renormalized matrix element and 
apply the corresponding counterterms obtained from the 



LO matrix element by inserting the relevant Z— 1 factors 
into the internal propagators and vertices. All the renor- 
malization constants we need are presented in Eq. (|2.5|) . 



We will get the renormalized vertex function , where 
(N) denotes the set of N external particles. The renor- 
malized matrix element is obtained from 



N 

no 



(2.6) 



(i) 

where Z R are the residues of the renormalized propaga- 
tors at the poles for all the particles under consideration. 
They are related to the residues of the unrenormalized 
propagators via 



y{i) _ y{l) 7 -\ 

^R ~ ^1 



(2.7) 



where the Zi are the respective external wave function 

renormalization constants. 

Working at the one-loop order, we note that in the on- 
to 

shell scheme Z R = 1. This is a direct consequence of 



the second condition in Eq. (|2.4[) . which effectively cuts 
off the external massive lines. For the case of external 
massless partons Zjf = 1. It is important to note that 
the gluon wave function renormalization constant Z3 is a 
mixture of two parts: the part which multiplies C £ (/j 2 ) is 
derived in the MS scheme, while the last term due to the 
heavy-quark loop is derived in the on-shell scheme. For 
this reason, this last term has to be omitted in Z3 when 
using it as an external field renormalization constant in 
Eq. (|2.7|) . Since in our case we have two gluon and two 
heavy-quark fields, we therefore obtain 



m r = r 



_ r ( JV ) ■ 



(2.8) 



The final result should not depend on which of the 
two ways has been chosen to do the renormalization. We 
have checked that, in both ways, one arrives at the same 
renormalized matrix element. 

In order to fix our normalization, we write down the 
differential cross section for gg — > QQ in terms of the 
squared amplitudes \M\ 2 . One has 



da 



1 <ps) 2 1 |M|2 



gg^QQ 2 s 4(1 - e) 2 d\ r " 'sg^QQ 



(2.9) 



where the n-dimensional two-body phase space is given 
by 



-2e 



(4tt) £ 



ins T(l -e) 



tu 



5(s+t+u)dtdu . 



(2.10) 

We explicitly exhibit the flux factor (4pip 2 ) 1 = (2s) 1 , 
and the spin (n— 2)~ 2 = (2 — 2e)~ 2 and color d^ 2 averag- 
ing factors for the initial gluons. Here — — 1 = 8 
is the dimension of the adjoint representation of the color 
group SU(A^c)- 
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FIG. 3: The t-channel one-loop graphs contributing to the gluon fusion amplitude. Loops with dotted lines represent the gluon, 
ghost, and light and heavy quarks. 



III. LEADING AND NEXT-TO-LEADING 
ORDER RESULTS 



At LO for gg — > QQ, we shall use a representation 
which differs from the one given in Refs. Hfl, EH. First 
note that there are only two independent color structures 
for this subprocess. The s-channel matrix element is a 
sum of two parts, each of which is proportional to one 
of the two independent color structures. We combine 
terms with the same color structures of the three (e.g. 
s, t, and u) production channels. Finally, we remove the 
heavy-antiquark momentum p± using energy-momentum 
conservation and use on-shell conditions for the gluons 
(Pi • ei = and P2 ■ e 2 = 0) and the heavy quark (u 3 fa = 
U3171). We then obtain the two color-linked LO matrix 



elements 

M LO ,t = iT b T a M/t, M L o,u = iT a T b M/u, (3.1) 
with 

sM = ^fcrfs + 2^p\t - 2~/ u p$t - 2YPz-9 - 2fag^t. 

(3.2) 

It can be verified that the function M is t «-> u symmetric, 
and consequently the color- linked Born amplitudes M L o,t 
and A/lo.u turn into one another under t «-» u. 

We then square the full Born matrix element M\_,o,t + 
A^lo.u and do the spin and color sums to obtain the LO 
amplitude, 
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FIG. 4: The s-channel one-loop graphs contributing to the gluon fusion amplitude. Loops with the dotted lines as in gl, h, jl, 
and j2 represent the gluon, ghost, and light and heavy quarks. The four-gluon coupling contribution appears in g2. 



where we have factored out a color-reduced Born term 
\M\ 2 , which reads 

\M\ = 8\ =— + 4 4— 

I s 1 s tu 

_e2(l-^)+ £ 2 } =B. (3.4) 

The expression in Eq. (|3.3|) for the LO amplitude 
agrees with the well-known result in n dimensions (see 
e.g. Ref. [I]). Note that, by using the prescription 
of Ref. [23|, we were able to avoid the introduction of 
ghost contributions which would otherwise arise from the 
square of the right-most three-gluon coupling amplitude 
in Fig. O In our case the prescription of Ref. [23| con- 
sists in the use of on-shell conditions for external gluons, 
i.e. pi • €\ = and p2 ■ £2 = 0, and the exclusion of the 
heavy- antiquark momentum via p^ = p\ +P2 — P3- When 
squaring amplitudes, we sum over the two helicities of 
the gluons using the Feynman gauge, i.e. we use 

e»{\y{\) = -g» v . (3.5) 

A=±l 

The use of the framework set up in Ref. [23J has the ad- 
vantage in the non-Abelian case that one can omit ghost 
contributions when squaring the amplitudes. Using the 
above on-shell conditions already at the amplitude level 
means that one takes full advantage of the gauge invari- 
ance of the problem when squaring the amplitudes. Thus, 
in general, the results for the different channels will not 
be identical to the ones which would be obtained using 
't Hooft-Feynman gauge throughout. 



Folding the one-loop matrix elements (see Figs. [3] and 
[4]) with the LO Born term (see Fig. [2]), one obtains the 
virtual part of the NLO result . 

As concerns the one-loop matrix elements, we shall use 
the one-loop matrix elements of Refs. [20L l2l| to com- 
pute the virtual NLO contribution up to 0(e 2 ) in terms 
of the coefficient functions (|T1D of the scalar master in- 
tegrals. However, in Ref. [20|, where expressions for the 
NLO matrix elements up to O(e ) are given, the val- 
ues for the scalar coefficient functions in terms of loga- 
rithms and dilogarithms are substituted directly. There- 
fore, we had to recalculate the corresponding expressions 
from Ref. pp| for the matrix elements in order to have 
a uniform result in terms of scalar coefficient functions. 
This has allowed us to retrieve and use relations between 
coefficients of the scalar coefficient functions in the result 
for different orders of the Laurent series expansion in e. 
We will comment on these relations later on. 

We also mention that we had to regroup and rear- 
range various terms in the one-loop amplitudes from 
Refs. [13, [2l| according to the three independent color 
structures in order to bring the pole terms into agreement 
with the form suggested in Ref. (24J. In the gluon fu- 
sion case treated here, there are three independent color 
structures in the one-loop amplitudes, e.g. x b T a ,T a T b , 
and 8 ab . As in the LO case, one also has to exclude the 
heavy-antiquark momentum P4 from the one-loop ampli- 
tude expressions. As a result of the above two steps, the 
pole terms of our new matrix elements became propor- 
tional to the LO color- linked amplitudes (|3.ip . In all our 
subsequent calculations, we shall use only these matrix 
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elements. 

The NLO virtual corrections to heavy-flavor hadropro- 
duction have been calculated before for the gg — > QQ 
case. Nevertheless, one cannot find explicit separate re- 
sults for the virtual corrections in the literature although 
Ref. provides analytic results for the combined "vir- 
tual+soft" contributions. We have therefore recalculated 
the virtual NLO contribution to gg fusion. In fact, we 
have calculated the virtual NLO results up to C(e 2 ). As 
it turns out, use of the expressions for the NLO virtual 
0(e 1 )- and O(e 2 )-contributions considerably simplify the 
presentation of the corresponding NNLO results in as 
much as they appear as important building blocks in the 
NNLO results. 

Next we fold the pole, finite, 0(e 1 ) and 0(e 2 ) terms 
of our NLO matrix element with the LO matrix ele- 
ment. In dimensional regularization, the trace evalua- 
tion in n = 4 — 2e dimensions will lead to terms of order 
0(e 1 ) and 0(e 2 ) when multiplied with the pole and fi- 
nite terms, as well as to the terms of 0(e 3 ) and C(e 4 ) 
when multiplied with the 0(e 1 ) and 0(e 2 ) terms of the 
squared amplitude, respectively. In the following we will 
disregard terms of 0(e 3 ) and 0(e 4 ) as they do not con- 
tribute to the finite part of the NNLO result. 

Before presenting our result for the NLO matrix ele- 
ment, we would like to comment on its color structure. 
We have decomposed our matrix elements according to 
the following three independent color structures: 



S ab Xr ( T a T fc) 

Tv(T b T a ) Ti(T b T a ) 



(3.6) 



2 

cLa . 
2 

Aa 

2 Ui 2 

At NLO, the final spin and color summed matrix ele- 
ment can be written as a sum of five terms: 



Tv{T b T a ) Tr(T°T h ) 



-C f 



N c 



)• 



\M\ 2 



Loop x Born 9. 



2 VC Re 



eW (1) (e)~ 



-W { ~ 

e 

E 2 W^ 



1) 



00 



(3.7) 



where C has been defined in Eq. (J273J). The notation 
| M | L 00 p x B orn means that one is retaining only the G(o? s ) 
part of \M\ 2 . 

The first two coefficient functions in Eq. 1|3.7|) have a 
rather simple structure: 

W^- 2 \e) =-4N c B, (3.8) 

W^\e) = d A B(^h + (C F ~ ^)(/ t + /„) 
u t 

+GF — ft + Cp — fu 

t u 

where B and B are the LO terms defined in Eqs. 
and (13.41). We have also introduced new functions, 



n 



i 



In- 



s m 



-hi — 

S 771 



2m 2 



2s(3 



In a;, 



f t = N c ln^+2N c ln- 



2C F - fa 



+ {2C F -N c )^r^ In x, 



s/3 



fu — ft\t* 



(3.9) 



where (3 — y/1 — 4m 2 /s is the heavy-quark velocity and 
(3 Q is defined after Eq. Q2.50 . 

One should keep in mind that the overall Born term 
factors B and B contain terms multiplied by e and 
e 2 . Therefore, if the expressions for B and B, given in 
Eqs. Q and (32]), are substituted in W^~ 2 ' and , 
we will obtain additional 0(s^ 1 ) and finite terms from 
the first two terms of Eq. (|3.7p . 

The third term in Eq. (|3.7[) reads 



^°)( £ )eF/ ' 



NLO ' 



(3.10) 



where we have constructed the following generic func- 
tions: 



NLO - W l 



0) 



with 



<I_A 

2 



= - 



(t <-> u)\ 
2BP 



(3.11) 



(3.12) 



The three functions Fi,F 2 , and F3 are defined as follows: 



(J) 



F 



with 



(a/ 



ea 



(<0 



e a 



F 



(i) 



jw _ {B2, B$, Ci, C 2 , C 2u , C3, Csu, C4, C5, Cq, 
D 1 ,D lu ,D 2 ,D 2u ,D 3 }W; 



1 



(3.13) 



with /«) = {1, B 2 , B B) d, C 4 , C 5) C 6 } 



0). 



F. 



U) _ 



EC 



(<0 



£ C 



with 



jU) — Bi,B 2 , B 5 ,Ci,C 2 , C3, C4, C5, Ce, 

For J = 1 one has = 1, otherwise 1^ = b[ 3 \ C [ 2 j) 
etc. In other words, the summation index I runs over 
the scalar integral coefficient functions, while the coeffi- 



cient functions aj,a T 



(<0 >*) 



etc. denote the explicit de- 



pendence on s,t and m 2 . These coefficient functions are 
presented in Appendix A. Note that index j takes the 
same value for all the coefficient functions in Eq. (|3.13[) 
as well as in similar equations that will follow. 

The additional subscript "u" in some of the scalar coef- 
ficient functions in the expression for f[ 3 ^ (such as C 2 ^ ) 
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is to be understood as an operational definition prescrib- 
ing a (t u) interchange in the argument of that func- 



tion, i.e. C. 



(o) 

2it 



a 

,U) 



(0) I 



etc. 



Note that W 2 is only contributed to by the renor- 
malization procedure. Of course, all the remaining 0(e) 
terms (e.g. W^\e) and W^ 2 '(e), as well as those coming 
from W^ 1 ) (e) / e and W^(e)) should be disregarded in 
the NLO final result in Eq. (|3.7|) . It is important to note 

that -PnloI e=0 is n °t formally the full finite part of the 
NLO result in dimensional regularization, but it results 
from folding the finite part of our original NLO matrix el- 
ement with the LO one. Another part of the finite result 
comes from the first two terms in Eq. (|3 . 7|) , as mentioned 
before Eq. (|3.10[) . However, one should realize that the 
first two terms in Eq. (|3.7p would be cancelled with the 
corresponding parts from the real bremsstrahlung dia- 
grams. Given the overall factor, Eq. (II. 2|) . the term F^£ Q 
evaluated for e = represents the finite part of the vir- 
tual one-loop NLO result. 

Our 0(e- 2 ), Oie- 1 ), and O(e ) NLO results in 
Eq. (|3.7[) were analytically compared with the corre- 
sponding results obtained in Ref. [lj, which were kindly 
provided to us in a Schoonschip format by the authors 
[25l | . We obtained complete agreement. 

The fourth term in Eq. (|3.7p is a result of folding the 
0(e) term of the matrix element with the Born term. Be- 
cause of the n-dimensional traces, one also obtains terms 
of 0(e 2 ) and C(e 3 ). As mentioned before, we will only 
retain terms of O(e) and 0(e 2 ). We have 



W {1) {e 



(i) 



NLO 



Ft 



(0) 



NLO.e' 



where 



F, 



(j) 



NLO,£ 



d A 



F 



-N C F 7 



{(C F + 
+ (*«-> u)} 



Nct_ 
2 s 



0) 



(3.14) 
(3.15) 



Here 



fP = £(4 £) + £ 4 s2) )/« 



with /W) = {B 2 , B 5 , C U C 2 , C 2u , C 3 , C 3u , C 4 ,C 5 , C 6 , 
D 1 ,D lu ,D 2 ,D 2u ,D 3 }^; 



^' ) = E(4 £) +4 £2) )^ ) I 

with /« = {1,B 2 ,B 5 ,C 5 } 
F« = £(^+ £ff fV (i \ 



0). 



(3.16) 



0). 



with JW) = {l,B u B 2 ,C 2 , C 5 ,C 6 ,-Di} y 

F 7 W =^(4 £) + £ 4 £2) )/«, 
i 

with /W = {1, B u B 2 ,B 5 ,d, C 2 ,C 3 , C 4 ,C 5 ,C 6 , 
D!,D 2 }W. 



The coefficients di,ex,gi,hj are presented in Ap- 
pendix B. Note that the first term in Eq. (|3.14p in nothing 
but the NLO term of Eq. (|3.10p with indices of the co- 
efficient functions of the scalar master integrals and the 
power of the logarithm that multiplies (3q shifted upwards 
by one. 

The last term in Eq. I|3.7p is a result of folding the 
0(e 2 ) term of the matrix element with the Born term. 
Because of the n-dimensional traces, one also obtains 
terms of 0(e 3 ) and 0(e 4 ), which are omitted as before. 
For the 0(e 2 ) terms we obtain 



W (2) () _ F (2) F (l) 



NLO.e 



F 



(0) 



NLO.e 2 



(3.17) 



where 



NLO,e 2 



= d A 



r< p0') 



+JV c F 1 ( f) + (t« M )} 



(3.18) 



Here 



i 

with jC0 = {&, C 2 , C 2u , C 3 , C 3u , Ca, Cfc, C 6 ,Di,D 1v , 
D 2 ,D 2u ,D 3 }W; 

1^ ={l,B 2 ,B n ,C 5 }^; 

^mfh^, (3.19) 
with JW) ={1,C 2 ,C 5 ,C 6 ,D 1 }«; 



F, 



with 



F 
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with /0) = {1, B 5 , Ci, C 2 , C 3 , C4, C5, C*6, -Di, Dg}^ ■ 

The coefficients kj,lj,mj,nj are presented in Ap- 
pendix C. We mention that the functions Fi,^, and Fg 
are (£ <-> u) symmetric. 



IV. SINGULARITY STRUCTURE OF THE 
NNLO SQUARED AMPLITUDE 

The NNLO final spin and color summed squared ma- 
trix element can be written down as a sum of five terms: 



\M\ 



Loop x Loop 



= Re 



-U(-*)( e ) + -LyH0( e ) ( 4.i) 



r y(- 2 )(e) + i^- 1 )( e ) + y(°)( £ ) 



where C has been defined in Eq. (|2~5)) . Note that Eq. P~T]) 
is not a Laurent series expansion in e since the coefficient 
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functions y( m ) (e) are functions of e as explicitly anno- 
tated in Eq. (|4.ip . It is nevertheless useful to write the 
NNLO one-loop squared result in the form of Eq. (|4.1|) 
in order to exhibit the explicit e structures. All five coef- 
ficient functions V^ m \e) are bilinear forms in the coeffi- 
cient functions that define the Laurent series expansion of 
the scalar master integrals (11.11) . Some of these coefficient 
functions are zero and some of them are just numbers or 
simple logarithms. In the latter case, we have substituted 
these numbers or logarithms for the coefficient functions 
y( m ) m the five terms above. This has been done for all 
the scalar coefficient functions that multiply poles, i.e. 
for scalar functions with negative subscripts /(~ 2 ) and 

as well as for the whole scalar functions B { 3 1 \ where we have introduced new functions 



The last term -27V c i^ L ) o in Eq. (|4~3| is obtained from 
folding the 0(e~ 2 ) singular term of the matrix clement 
with its finite part, while the remaining parts result from 
folding the single poles. Note that when one substitutes 
the Laurent expansions for B and F^q, one gets addi- 
tional 1/e poles and finite terms in Eq. (|4.3|) . 

The structure of the fourth term in Eq. (|4. 1[) is some- 
what more complicated. One has 

V^(e) = ^ln(4)V^He)+s[ a) - 2N c W^(e), 

(4.5) 



and Bf. 

We found that a significant part of the NNLO results 
can be expressed in terms of the e expansion of the NLO 
contribution. In particular, we will need the NLO expan- 
sion up to e 2 . Therefore, in this section, we will make full 
use of the results derived in Sec. IIII1 

Before proceeding further, we note that there are no 
additional color structures appearing in the NNLO cal- 
culation for gg fusion in addition to the ones already 
presented in Eq. (|3.6p : they are just linear combinations 
of the ones in the NLO case. This is in contrast to the 
qq subprocess, where the NNLO color structures exhibit 
much higher complexity and richness [l5j relative to the 
NLO ones. 

The two most singular terms in Eq. (|4.1[) are propor- 
tional to the Born B and color-reduced Born B terms 
defined in Eqs. (|3.3[) and (|3.4p . respectively. One has 



V { - 4) {e) = AN'^B, 



V { - 3) (e) = -2N c W { -' } {e), 



(4.2) 



where W i --~ 1 \e) is given in Eq. (|3.8p and is nothing but 
the full coefficient of the single-pole NLO result. 
For the 1/e 2 term we obtain 



V^He) = d A B S -\f s \ 2 + L Cl 



u 

r s (- t ft + ^fu) + (c F -^)f;f, 



(4.3) 



where the functions fsifti an d f u above are the same 
as those in Eq. (|3.9p . but now with the imaginary parts 
retained, i.e. one has the following replacements: 



In — - — > In — - — in, 



lnx — > lnx + itt. 



(4.4) 



This reflects the fact that, contrary to the NLO calcula- 
tion, one has to keep the imaginary parts in the NNLO 
calculation as emphasized in the Introduction. It should 
be clear that the completion (|4.4p has to be done every- 
where in the NNLO calculation whenever the logarithms 
(14. 4p appear in bilinear forms multiplying complex func- 
tions. 



= (HF® + LtF^ + L^fP + (t ~ „)) , 



with 



Li = 2f s f t f u , 

s s 

u t 

L 2 = 2/ 4 - 2C F -f t - {2C F - Nc)~f u 
s ' s 



(4.6) 



(4.7) 



The first two terms in Eq. (|4.5p arise from folding the 
single-pole terms in the original matrix element with its 
finite O(e ) part. The last term is due to the interference 
of 0{e~ 2 ) x 0{e) terms in the original matrix element. 
This pole term is due to the Laurent series expansion 
of the original matrix element and cannot be deduced 
from the knowledge of the NLO terms alone. The func- 
tion IfWfe) is defined in Eq. (|3.14p . while the functions 
F^\F^ j \ and FP are given by Eq. {333]). 

When one substitutes the Laurent expansions for^ (0) , 

i*2°\ ^3 j an d W^\e), one gets finite and 0(e) terms 
in Eq. (|4.5[) . However, since we are only interested in 
the Laurent series expansion up to the finite term, these 
0(e) contributions can be omitted as before. 



V. STRUCTURE OF THE FINITE PART 

In this section, we present the finite part of our result. 
In the course of our calculation, we have made full use of 
the results presented in Sec lIHi e.g. of our detailed study 
of the NLO structure of the Laurent series expansion up 
to 0(e 2 ). As a consequence, we can present a large part of 
our results for the finite part in a surprisingly concise and 
closed form. We decompose the finite part into several 
pieces, as 



y (0) (e) = Re V{[ 



r(0) 



V, 



(o) 



(0) 



(5.1) 



The first two terms originate from the interference of 
the ©(e" 1 ) x 0(e) and 0(e~ 2 ) x 0(e 2 ) pieces of the ini- 
tial matrix element, respectively. Each of them can be 
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conveniently presented in a very compact form: 



+s>» +sf, 



-c. 



where we have introduced one more function, 



(5.2) 



5. 



(.;) 



0} 



(5.3) 



Similarly, for the second term in Eq. (|5.ip . we write 



y(0) 



~2N c W {2) {e), 



(5.4) 



with W {2) (e) defined in Eq. (pTT7|) . Note again that the 
0(e) and C(e 2 ) terms in the above expressions for 

and V22 can t> e disregarded. We mention that the scalar 
coefficient functions with the superscript "2" above in- 
volve multiple polylogarithms of weight and depth 4. 

We emphasize that the quasifactorized forms of all the 
expressions given in this paper hold only when one retains 
the full e dependence in the Born and NLO terms. 

The last term in Eq. (|5.1[) comes from the square of the 
O(e ) term of the matrix element, which can be written 
as 



y(0) 
v 00 



Tft 

-A)ln(— ) 



7 (0) 
NLO 



1. 



(0) 



Y . 



(5.5) 



where F^ Q and W^ 0) are given in Eqs. |3~TT|) and 1(532)) . 
We found that the last term Y in Eq. (|5.5[) also possesses 
the quasifactorization properties discovered in a recent 
paper [l5| . For instance, the result can also be written 
down as a sum of bilinear products, where each of the fac- 
tors are linear combinations of scalar integral coefficient 
functions multiplied by some combinations of kinematic 
variables. However, because of the great number of Lau- 
rent structures appearing in the original matrix element 
for the gg fusion subprocess, the length of the final ex- 
pressions does not allow us to present the results in this 
paper. Also, we were not able to find the optimal way to 
organize the different contributions in Y as in Ref. 15], 
as not all the powers of common numerators and denom- 
inators cancel out. Therefore, we have opted to supply 
the results on the finite term Y in a separate electronic 
file. 

In the finite contribution of Eq. (|5.ip , one notices the 
interplay of the product of powers of e resulting from 
the Laurent series expansion of the scalar integrals [cf. 
Eq. (jl.ip ] on the one hand and powers of e resulting from 
doing the spin algebra in dimensional regularization on 

the other hand. For example, for the finite part one has 

(— 1) (0)* 

a contribution from B] as well as a contribution 



from C^B^* 



Terms of the type Cg ^ Bj * , where 



the superscripts corresponding to e powers do not com- 
pensate, would be absent in regularization schemes where 
traces are effectively taken in four dimensions, i.e. in the 
so-called four-dimensional schemes or in dimensional re- 
duction (DRED). 

We emphasize that all our factorized results given in 
this paper [except for the expression for Y in Eq. (|5.5[) ] 
take up about 22 Kb of hard disk space. This has to 
be compared with the length of the original, untreated 
FORM output. The original computer output for the 
corresponding one- loop squared cross section of the gg — > 
QQ subprocess turned out to be very long and took up 
about 85 MB of hard disk space. Therefore, the reduction 
is of the order of 10 3 -10 4 in the present case. 

As a final remark we want to emphasize that we have 
done two independent calculations using REDUCE (2(| 
and FORM 1271 when squaring the one-loop amplitudes. 
The results of both calculations agree. Casting the re- 
sults into the compact forms presented in this paper was 
done with the help of the REDUCE Computer Algebra 
System. 



VI. 



CONCLUSIONS 



We have presented analytical 0(a A s ) NNLO results for 
the one-loop squared contributions to heavy-quark pair 
production in the gluon-gluon fusion reaction. The cor- 
responding result for photon-photon fusion has already 
been presented in Ref. [l7|, while results for the photon- 
gluon fusion process can be obtained from Ref. [21| after 
some color factor adjustments. As concerns hadropro- 
duction of heavy quarks, the results of the present paper, 
together with a recent publication on qq production [1 51 ] , 
complete the derivation of the one-loop squared contribu- 
tions to the hadroproduction of heavy quarks at NNLO 
with the heavy-quark mass dependence fully retained. 
Our results form part of the NNLO description of heavy- 
quark pair production relevant for the NNLO analysis of 
ongoing experiments at the TEVATRON and the LHC. 

A large part of our analytical results are presented in 
a very compact form. The singular contributions propor- 
tional to £~ 4 ,e~ 3 , and e~ 2 are entirely given in terms of 
LO and NLO contributions, whereas the e _1 contribu- 
tions contain some true NNLO structure in addition to 
LO and NLO structures. Since the LO and NLO terms 
are themselves expanded in Laurent series, this implies 
that our singular contributions are not true (in a math- 
ematical sense) Laurent series in e. We believe that our 
representation of the singular contributions has struc- 
tural advantages in as much as it will be simpler to match 
our singular structures onto the singular structures of the 
other classes of contributions. Also, our representation 
is convenient if one wants to convert our expressions to 
different regularization schemes such as DRED (see e.g. 
Ref. [28(]). If needed, our singular contributions can eas- 
ily be converted into true Laurent series expansions since 
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our expressions are very compact. 

Because of our representation of the singular parts, we 
obtained quasifactorized expressions for a large part of 
the finite contributions. Writing our analytical results in 
factorized forms led to a reduction of the length of the 
original output by a factor of 10 3 -10 4 , which will lead to a 
dramatic reduction of the CPU time needed in numerical 
evaluations. 

The present paper deals with unpolarized gluons in the 
initial state and unpolarized heavy quarks in the final 
state. Since our results for the original matrix elements 
contain the full spin information of the process, an exten- 
sion to the polarized case with polarization in the initial 
state and/or in the final state including spin correlations 
would be possible. 

Analytical results in electronic format for the coeffi- 
cients given in the Appendices as well as for the term Y 
in Eq. ()5.5j) are readily available p9l |. 
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APPENDIX A 

First, we write down a few abbreviations that we use 
throughout the paper: 

ff = yl — 4m 2 /s, D — m 2 s — tu, 

z 2 = s + 2t, z 2u = s + 2u, (Al) 
Zf = 2m 2 + t, z u = 2m 2 + u. 

Note that D in Eq. (|A1|) is not the space-time dimension. 

Here we present the expressions for all the coefficients 
or, & r, ci appearing in Eq. (|3 . 13[) : 



a B2 = 16D/{sff 2 ) , 

as 5 = — ob 2 , 

a Cl = 4(8m 4 - z 2 /s{2m 2 - s + 2m 2 j ff 2 )) , 

a C2 = 8t/s(4m 2 z t + 2st + t 2 ) , 

a c 2u = a c 2 ( f <-> u ) i 

a Cs = 8t/s{Am 2 z t +tz 2 ) , 
a c 3u = a Ca (t<->u), 

a Ci = 4(4to 2 s + 3s 2 - 8tu) , 

a Cs = 4(8m 4 - 3s 2 + 2tu) , (A2) 
a Ce = -4f3 2 (2m 2 s + s 2 + 2tu) , 
a Dl = 4(2to 2 (2L> + sztP 2 - t 2 [3 2 ) + s 2 t0 2 + t 3 ) , 
a D lu = a Dl (t <-> it) , 

a D2 = 4(8to 2 L> - stuff 2 + 2t 2 /s(t 2 + u 2 )) , 
a D 2u = a D2 (t <-> u) , 

a D:i = 8(8rn 2 D - 8mHu/s - stuff 2 + 2t 2 u 2 /s) ; 

ag = 4(2 s -z 2 /( s/ 3 2 )), 
0) _ _ 0) 

a B 5 ~ a B 2 ' 

ag = 2(ff 2 (s 3 (8m 2 + s)-8t 2 u 2 )/D 
+ 16m 2 /s(s 2 +tu + D/ff 2 ) , 

aP 2 = -4t 2 (10 - t/s(2tuff 2 + 2s 2 - M 2 )/D) , 

a cL = a c 2 ) (* ^ u ) - 

= -4t 2 (6 - t/s(2tuff 2 + 2s 2 - Ast - 5t 2 )/D) , 

a cL = ^cli 1 ~ u ) ' 

a { cl = -2(s 3 (2m 2 - s) - 8tu(m 2 s + t 2 + u 2 ))/D , 

ag = 2(/3 2 (s 3 (6m 2 - s) + U 2 u 2 )/D (A3) 
+ 8s 2 - 12m 4 z|/L») , 

a g = 2(ff 2 (s 3 (8m 2 ~s)~4t 2 u 2 )/D 
+ 8m 2 s~ Am^zl/D) , 

a%\ = -2t(2s 2 ff 2 ~ s 2 tff 2 (2m 2 zl/s 2 + 4m 2 + t)/D 
+ 2tz t + 2s 2 ) , 

a r>l u = u ) . 
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a% ] = -2t(2s(s -u) + t 2 (s 2 + 8tu - 8u 3 /s)/D) , 



,00 



= Atu{As - t 2 u 2 /s 2 (8m 2 - 7 s)/ D) ; 



> 2 



l C 3 

■SI 

> 2 



(e 2 

in 

L>1 



a Dr = 



= o, 

- o, 

= 8s(s- 2m 2 zl/D) , 

= -8t 2 (3u/s + t(m 2 - u)/D) , 

= a { S 2 \t <-» u) , 

= 8t 2 {2 + tu(l-3t/s)/D) , 

= a c 3 } ( f ^ u ) ' 



-16s 2 tu/^ . 
-4s(s + 2m 2 z 2 /D) , 
As(s-2m 2 zl/D) , 

-4st(2m 2 - s + t(fi 2 tu + m 2 z 2 /s)/D) , 

4st(s- U 2 u/D) , 

a^ 2 \t <-> u) , 

-8tw(s + 3tV/(sL>)) ; 



(A4) 



6i = -16/3z 2 /s(m 2 (n ; + 1) + (2C F - N c )3D/(s/3 2 ) 
- N c (m 2 + L>6(10m 2 - s)/(s 2 [3 4 ))) , 

b B2 = -8z 2 /s 2 (8m i - (2C F - N C )D(2 + 1//3 2 )) , 

&b 5 = -^c8z 2 (£»(16TO 2 -s)/( s /3 4 )+tw)/s 2 , 

6 Cl = -^ c 16m 2 J Dz 2 (8m 2 + s)/( S 3 /3 4 ) , 

6 C4 = N c 4z2(D-2tu)/8 , (A5) 

&C 5 = -32m 4 2: 2 /s , 

b Ce = -{2C F - N c )16Dz 2 (2m 2 - s)/s 2 ; 

&< e) = m/3z 2 (tu(ni + l) + {2C F - N C )3D/P 2 

- N c (36m 2 D/(sf3 4 ) - tu(4m 2 - 7s)/(sf3 2 )))/s 2 , 



b% = 



$ = 
$ = 



cs 2 

CJ3 5 

cc 2 



c c 3 

C C 6 



8z 2 (8m 2 tw/ s + (2C F - N c )(2tu - D/j3 2 ))/s 2 , 
N c 8z 2 (3m 2 z 2 /(sl3 4 ) - 2{D + 2m 2 tu/ s) / (3 2 ) / s 2 , 
N c 16m 2 z 2 (3D/[3' i + 2tu/(3 2 )/s 2 , 
N c l2tuz 2 /s , (A6) 
32m 2 iiiz 2 /s 2 , 

-{2C F - N c )mtuz 2 (2m 2 - s)/s 2 ; 



ci = 16(C F (L>/3 2 (8m 2 T/t 2 + 2) - D(6z t /t - 2 - t/s) 
+ 2m 2 {Az t {m 2 /s - 1) - m 2 ) - D(l + 4i/s)//3 2 ) 
- N c {D2m 2 (2D + tu)/(st 2 ) - 2m 2 tu/s 
-D4m 2 {s + U)/{s 2 f3 2 )))/T , 

c Bl = 16{C F (2m 2 f3 2 (T-2s-D(2T + t)/t 2 ) 
+ D(3z t /t + t/s) - 2m 2 u(2 + ht/s))/T 
+ N c 2D(D/s-t)/t 2 ) , 

(2C F - N c )WD/(sf3 2 ) , (A7) 
N c 8(-8m 2 D/(s 2 f3 2 ) - t(3 2 + t 2 z 2 /s 2 ) , 



Af c 8(t 3 + u 3 - 4t 2 T - sD/0 2 - s 2 (3 2 {m 2 - t))/t 



-{2C F - N c )l6(2m 2 z 2 (m 2 s/t - z t ) 
+ t(s 2 +t 2 ))/s , 

N c l6t(AD/s-t/3 2 +s) , 

A c 4(-s 2 /3 2 + 3z 2 (m 2 s - t 2 )/s - 3su + 2t 2 ) , 

(2C F - Nc)8(2T(2m 2 + s) - u 2 ) , 

-{2C F - N c )8(Am 2 D/s - 4m 2 t/3 2 + 3tz t - z 2 2 ) , 

c Dl = -(2C F - N c )8{m 2 s 2 (3 A - 2m 2 t(3 2 {s - t) + st 2 (3 2 
-t 3 - sD) , 

Cd 2 = N c 8(8m 2 D - stu/3 2 + 2t 2 (t 2 + u 2 )/s) ; 

c[ E) = 16(C F (D(16m 2 D/{st 2 ) - 24m 4 /t 2 +4- t/s 
+ 2t/{s0 2 )) + 2m 2 (4m 2 - 6* - 9t 2 /s) 
+ 4m 2 t 2 z 2 /{s 2 p 2 ))/T 
+ N c 2{2m i s/t 2 + t + Dz 2 /{st) 

- D(4m 2 + 3s)/{s 2 [3 2 ) + tz 2 /(s[3 2 ))) , 

eg = 16{C F (4m 4 D/t 2 - QTD/t-2m 2 D/s -tD/s 

- bm 2 z t +t 2 ) 

- A c (2m 2 J D 2 /(st 2 ) + 2W/s - m 4 + t 2 ))/T , 



-(2C F -N c )8t(2 + z 2 /(s/3 2 )) , 



(A8) 
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N c 8(2D/(3 2 + 6m 2 z 2 //3 2 - 3t 2 - 2t 3 /s)/s , 

-N c A(2m 2 z 2 /s - 'is 2 - At 2 - 4m 2 (4z t 
+ tz 2 /s)/(3 2 + 2tuz t (As + 3i 2 /s + u 2 /s)/D 
+ st 2 (2t/3 2 - z 2 )/D) , 

(2C F - N c )8{f3 2 t(6sD - Am 2 tu - st 2 ) 
+ 2D(2m i s/t - sz t + Dt/s - t 2 ) 

- Am 2 t 3 z 2 /s)/D , 

-N c 8t 2 {As/t + 14 - st(Af3 2 - 8tT/s 2 + 5)/£>), 

~N c A(2s 2 + 2t 3 /s - 2su + m 2 st(9s + It) /D 
-t i (9 + 8t/s)/D) , 

-{2C F - N c )A(2m 2 z 2 /s - 2t 2 

- Asf3 2 (s - m 2 tu/D) + t 2 (8m 2 t + s 2 )/D) , 

-{2C F - N c )Ast(7f3 2 + 2u 2 f3 2 /{st) + 2s /t + 5 
+ /3 2 (2m 2 z 2 /s-3st-4t 2 )/L>) , 

-(2C F - N c )Ast 2 (2sf3 2 /t + 2s/t + 2t/s 
+ 4m 4 zi/(s 2 D) 

- P 2 (6m 2 s - Am 2 tu/s + st)/D), 
~N c At(As 2 + 2st + t 2 (z 2 + \2tu - 8u 3 /s) / D) ; 



-16(C F (2m 2 (6D/t + Am 2 + t) 
-Am 2 (3D + tz 2 )/(s(3 2 ) + D/f3 2 )/T 
-N c 2(u-Am 2 z u /{sf3 2 ))) , 

16(C F z t (3D/t + 2m 2 )/T- N c 2m 2 ) , 

-(2C F - N c )32m 2 z 2 /(s(3 2 ) , 

-N c 6Am 2 z 2 /{sf3 2 ) , 

~N c 8s{At + 2m 2 t{s + At + z 2 /(sf3 2 ))/D 
+ s/f3 2 ) , 

(2C F - N c )lQ{2m 2 s + tu + 2m 2 t(m 2 s 
-t 2 )/D), 

N c 16t{s + 2t(m 2 s + tu)/D) , 
N c 8s(s + 2t(m 2 s + tu) / D) , 
(2C F - N c )8s(u - 2m 2 tz 2 /D) , 
-(2C F - N c )8(2m 2 z 2 (m 2 s - t 2 )/D + su) , 



(A9) 



= -(2C F - N c )8t(sz u + 2m 2 tz 2 /D) , 
c^J = N c 8st(s - At 2 u/D) . 

APPENDIX B 

In this Appendix, we present the expressions for all the 
coefficients di,ex,gi, hi appearing in Eq. (|3 . 1 6[) : 



41 
41 
41 

41 
41 
41 
41 
41 
41 
41 

41 
41 
41 
4 £) 

L>3 

with 



d 



B-, 



Br, 



41 



2s(Am 2 + z 2 /(sf3 2 ))/(tu) , 

-41 , 

s{s 2 (3 2 {8m 2 s + s 2 + 2D) + As 2 D 
- 16m 2 D 2 /(s(3 2 ) - 8m 4 z 2 )/(tuD) , 

2t(2u(D + m 2 s) + st 2 + n d t/s)/(uD) , 



d<g(t «->«), 

-2t(2m 2 s 2 + st 2 - K d t/s)/(uD) , 
<©*<-«). 

s 2 (n d + 3sD - s 2 (m 2 - s))/(tuD) , 
-s(k c + 2m 2 sz 2 )/(tuD) , 

— SK c /(tuD) , 

st(z 2 + [3 2 (n d - s 2 (m 2 - t))/D)/u 
dg(t «->«), 

st(n d + sD + s 2 (m 2 - t))/(uD) , 

dg(t<-u), 

2tun d /{sD) , 

4L» 2 - s(3 2 {8m 2 s 2 - 8m 2 tu - s 3 ) , 
10m 2 s 2 — 8m 2 tu — 3stu ; 



-8m 2 z 2 / (stu(i 2 



= -d 



B 2 



(Bl) 



-s((22mV - I6m 2 tu + s 3 )s/3 2 + Am 2 sD 
- 16m 2 D 2 /(sf3 2 ))/(tuD) , 

2t(6sD + Am 2 sz 2 + t 2 z 2 - n d t/s)/{uD) , 
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4 2 
4 2 

( £ 2 

-Di 
(e 2 

(e 2 
(e 2 

(e 2 
"D 3 

with 

K c 



dg } (* ^ u ) . 

-2t(K d t/s - s(2m 2 s - 2si - t 2 ))/(uD) , 
dg, 2) (*<->«), (B2) 
-s 2 ( Kd + sD + s 2 (m 2 + s))/(tuD) , 
-s(k c - 2to 2 s(4L> + z 2 ))/(tuD) , 

— SK c /(tuD) , 

st(2uL> - /3 2 (K d - 4siu - si 2 ))/(uL>) , 

-st(K d + s(2m 2 s - 2st - t 2 ))/(uD) , 
d%^{t <-> u) , 
-2tuK d /(sD) , 

4iwL> + s/3 2 (18m 2 s 2 - 16m 2 tu - s 3 ) , 
18m 2 s 2 — 16m 2 tu + stu ; 



e B2 = 3 (8m 2 + s)kik 2 , 
eg = 18m 2 sKiK 2 , 
4? =41/^, (B3) 

= 41 / K 2 , 



= 2s(n; + 1)kiK2 

4V = 41 / k 2 , 

with 

ki = 8z 2 /(9s 2 ) , 
n 2 = —m 2 s/{tu) 



Next, wc introduce common factors that appear in the 
various coefficients gi and hj. They are multiplied by 
one power of e and read 



s b2 = 2m 2 s - tz 2 /(s(3 2 ) 



(B4) 



s C 2 = ts C 5 — 4m suD , 

s c5 = 2D(D + s(8m 2 +t)) + 2st0 2 {2m 2 u - t 2 ) + st 2 z 2 , 
s c6 = D(3s/3 2 + z 2 ) + /3 2 (6m 2 s 2 - 8m 2 tu + s 2 t) . 

For the coefficients gg , we have 

g[ e) = 8{2m 2 s0 2 {4sT 2 /t + tz t - 2tu) + D(10m 2 u - bsz t 

- 2t 2 ) - 2m 2 t(s 2 + u 2 ) - D 2 m/(sf3 2 ) 

- 3Dt 2 z 2 /( S p 2 ))/(t 2 uT) , 



<?g = -8(D(4m 2 u - st) - t 3 (2sf3 2 + 3z t ))/(t 2 uT) , 

41 = ~$s b2 /(tu) , 5 g = 8s c2 /(Dtu) , (B5) 

<?g = 4ss c5 /(Dto), 5 g = 4 Mc6 /P«)> 5 g = -fcg. 

Finally, we introduce factors that are common to various 
coefficients gi and hi that are multiplied by two powers 
of e: 

c b2 = 2m 2 u + D , (B6) 
c c2 = 2D(4m 2 su - tD - st(l7m 2 + 3t)) 

- 4st 2 /3 2 (2m 2 u - t 2 ) + st 2 (3sz t + 4m 2 z 2 ) , 
Cc5 = 2D(20m 2 s - st + t 2 ) + 2st/3 2 (4m 2 u + st - 2t 2 ) 

+ 5st 2 z 2 , 

c c6 = 2D{s[3 2 -u)+ /3 2 (16m 2 (.s 2 - tu) + s 2 t) . 



(el 

For the coefficients g\ , we have 



9i 



(e ] = 8{D 2 \2t/{s0 2 ) + D 2 16m 2 /t + D4(2sz t - tu) 



- Dt(Um 2 + M)//3 2 - 2m 2 {l2m 2 s 2 T/t + 5t 3 ) 

- 4m 2 t 2 z 2 //3 2 )/(t 2 uT) , 



g% } - 8(L»(4^ t /t 2 + l/u) + 2m 2 (t/u - 2))/T , 



9b 



(e } -16z 2 c b2 /(stu/3 2 ) , g<£ 2 ] = 8c c2 /(Dtu) 



9 ( c 5 ) = -4sc c5 /(Dtu) , g£ 6 ' = -4sc c6 /(Du) 

9D! - l 9c e ■ 

For the remaining coefficients hi, we get 



Je 2 ) _ 



(B7) 



h[ e) = 8(tz 2 (2m 2 + D/(sf3 2 ))/f3 2 ~ m 2 {4sD /t - 2s 2 
+ tz 2 /9))/(t 2 u) , 



-8(3m 2 s + te 2 )/(» 



41 = ^ b2 /{tu) , 



(B8) 



h%\ = 4(2D(8m 2 t + s 2 )/(s 2 /3 4 ) + 4/3m 2 {s - u) 

- sz t /p 2 )/{tu) , 

hQ = -2s(s 2 f3 4 t/D + 2s(3 2 {2m 2 z 2 -tu)/D 

- 8m 2 s 2 z t /(W) - D6/t + 8t + 2tz 2 /{sf3 2 ) 

- DAz t /{stp 4 ))/u , 



41 = -^ c2 /{Dtu) , = 2t/sh { Ci , 

/ig = 2s{8m 2 u 2 /D + 4s - st(8m 2 + s)/D)/u , 
/ig = -2ss c5 /(Dtu) , /ig = -2 SSc6 /(^ U ) 



h {e) 



-th 



(e) 
C 6 ' 



41 
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(e 2 ) 
Bi 

B 2 

(e 2 ) 
B 5 



8(4m i s 2 /t - s 2 z t - W/9st 2 - t 3 /3 - 2/9t 4 /s 

- Dtz 2 j{s^) + t{uz 2 + 8m 2 /sD)/[3 2 )/{t 2 u) , 

16(Dz 2 +t 2 u)/(t 2 u) , 

8z 2 c b2 /(stuf3 2 ) , (B9) 

-8z 2 (4m 2 D/(s 2 f3 4 ) + tu/(3s) 

- 2m 2 u/{s[3 2 ))/{tu) , 

-2s{AQm 2 s/t + (8m 2 suf3 2 + 20m 2 su + I6m 2 t 2 

,2„2 , 



sH)/D + 2(2mV/f + Am 2 t - s 2 )/{sf3 2 ) 



+ 8m 2 D/(stf3 4 ) + 4m 2 (l/s 2 + /3 2 /D)z 2 //3 4 )/u, 



= -II,. 



-4c c2 /(Dtu) , hff = 2t/shff , 
-2s(20mV + AstuP 2 + stz 2 )/{Du) , 
2sc c5 /(Dtu) , h^J = 2sc c6 /{Du) 



(z 2 ) h (z 2 ) _ _ th ie 2 ) 



APPENDIX C 



In this Appendix, we present the expressions for all the 
coefficients fc/, To/, rij appearing in Eq. (|3.19[) using 
the following abbreviations: 

Kcl = -s(3 2 (18m 2 s 2 - s 3 + 2(2m 2 + s)zj) - 8m 2 sD , 
Kc6 = -sl3 2 zl + 16sD + 6stu . (CI) 

We have 



k {£2) - 

k {£2) - 

K c 2 ~ 

k {£2) - 

K c 3 - 



k {£) 

K c 3u 



— SK c i/(tuD) , 

2t(n c6 t/s + 2m 2 su - s 2 z t + t 3 )/(uD) 
2t 2 (K c6 /s + 2s 2 + u 2 )/(uD) , 



(C2) 



fc^ 2) = s 2 (k c6 + 2s(u 2 - st))/(tuD) 



Se 2 ) _ 



-s(k c1 - 2s 2 (2D - tuft 2 - zl))/(tuD) 



k c? = s 2 2 K c6 /{tuD) 



& 2 



A? 2 

"D 2 
'D, 



stf3 2 (n c6 - suz 2 )/(uD) , 
st(n cQ + s 3 - s 2 t)/(uD) , 
2tu(n c6 /s - st + u 2 )/D ; 



lf ] = 4/3s(n; + 1)kiK2 , Ig? = (^m 2 + 5s) Ki k 2 , 



lg ^ = 5sniKlK 2 , 



,(e 2 ) 



18m sk\k 2 



(C3) 



= 32(2m 2 (2mV/£ 2 - s 2 /t + 2t- 8DT/t 2 
-s/3 2 (2D/t + m 2 ) + D4U/(s/3 2 ))/(tu) , 



l c 2 



= 2z t /( S/ 3 2 )mg ) , 
z t /(tp 2 )m ( £ , 
4s(3 2 (k c6 /u- sz 2 )/D , 



(C4) 



= — tm 



(e 2 ) . 



> 2 ) 



,(- 2 ) 



M 2 ) 



= -16m 2 (18s 2 z t /t 2 + 82/3s + 2/M - 9sz 2 /(s{3 2 ) 
+lUz u D/(s 2 p A ))/(9tu) , 

= I6m 2 z 2 /(9tu) , 



= 2z 4 /( S /3 2 )n^ e 



^ 2 ) 



2s(k c6 + s 3 - s 2 t)/{Du) , 
— 2s0 2 (h c q — suz 2 ) / (Du) 



(C5) 



= —in 



(e 2 ) 
c 6 



M 2 ) - 



-tn 



(e 2 ) 
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